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In this paper, we derive and study two versions of the short pulse equation (SPE) in (2 + 1) -dimensions. 
Using Maxwell's equations as a starting point, and suitable Kramers-Kronig formulas for the permittivity and 
permeability of the medium, which are relevant, e.g., to left-handed metamaterials and dielectric slab waveg- 
uides, we employ a multiple scales technique to obtain the relevant models. General properties of the resulting 
(2 + l)-dimensional SPEs, including fundamental conservation laws, as well as the Lagrangian and Hamilto- 
nian structure and numerical simulations for one- and two-dimensional initial data, are presented. Ultrashort 
ID breathers appear to be fairly robust, while rather general two-dimensional localized initial conditions are 
transformed into quasi-one-dimensional dispersing waveforms. 
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I. INTRODUCTION 

Ultrashort pulses, having a duration of a few optical cy- 
cles, have been the subject of intense study over the last years; 
this is due to the fact that they find many applications in 
various contexts, ranging from light-matter interactions, har- 
monic generation, attosecond physics, nonlinear optics, and 
others 1 1]. A theme of particular interest related to ultrashort 
pulses, is their evolution in nonlinear media characterized by 
an intensity-dependent refractive index. In that respect, mod- 
els describing ultrashort pulses in nonlinear media, as well as 
their systematic study, have attracted much attention; see, e.g., 
Refs. I042]- Prominent examples are the modified Korteweg- 
de Vries (mKdV) equation yO, the sine-Gordon (sG) equation 
H, 0], combined mKdV-sG equations among others. 

Note that most of these works refer to the one-dimensional 
(ID) setting; in the two-dimensional (2D) one, pertinent stud- 
ies are related to few-cycle solitons described by the general- 
ized Kadomtsev-Petviashvilli (KP) equation llOTl and collapse 
dynamics of ultrashort spatiotemporal pulses 111 ill . 

Another model describing ultrashort pulse dynamics is the 
so-called short-pulse equation (SPE), which was first derived 
in the context of nonlinear fiber optics lH2h and later in the 
context of nonlinear metamaterials lH3ll . From the physical 
point of view, the interest in the SPE model arises from the 
fact that its few-cycle pulse solutions have been shown to 
compare more favorably to the ones of the original Maxwell's 
equations, as compared to pertinent solutions of the more tra- 
ditional nonlinear Schrodinger (NLS) model lfl2i [l4ll . Fur- 
thermore, this model is also interesting from a mathematical 
point of view, due to the existence of an infinite hierarchy of 
conserved quantities lH5ll . its connection to the sG model and, 
thus, to its complete integrability lfl6ll . The SPE admits var- 
ious types of solutions, including singular soliton solutions - 
the so called loop solitons II 1711 - as well as other non- singular 
solutions, such as peakons, breather- and periodic-type wave- 
forms (131 [l7l - [l9ll . Note that, recently, wave-breaking phe- 
nomena 120], as well as the global well-posedness question 
ll2lll of the SPE were also investigated. The above volume of 
work refers to the ID setting; to the best of our knowledge, 



the SPE has not been considered or analyzed so far in 2D. 

In this work, we derive and study two different versions of 
the SPE in (2 + l)-dimensions, namely the SPE-I and SPE- 
II. In particular, starting from Maxwell's equations, and as- 
suming general Kramers-Kronig or Sellmeier formulas for the 
permittivity and permeability (see, e.g., Ref. 112211 ), we use a 
multi- scale expansion method to obtain the SPE-I and SPE-II 
models. Then, we study the general properties of each model, 
present the Hamiltonian, Lagrangian and momenta, and also 
obtain zero-mass constraints that are used in the simulations 
(and, specifically, in the preparation of the initial data). Next, 
we explore the dynamics of ultrashort pulses in 2D, employing 
as initial conditions, either the breather solution of the under- 
lying ID SPE or a waveform localized in 2D; our purpose is 
to investigate if ultrashort pulses are prone to transverse insta- 
bilities (induced by the presence of diffraction in our models) 
and also identify purely 2D structures that can be supported 
by the SPE-I and SPE-II. We find that the ID breathers are 
stable in the 2D setting, while 2D initial data are gradually 
transformed into quasi- ID waveforms, reminiscent of the ID 
solutions. 

Our presentation is structured as follows. In Sec. II, we 
derive the SPE-I and SPE-II models. Sections III and IV are 
devoted to the general properties and numerical study of SPE- 
I and SPE-II, respectively. Finally, in Sec. V, we summarize 
and discuss our conclusions. 



II. DERIVATION OF 2D SHORT-PULSE EQUATIONS 

We consider the propagation of TE Z (transverse-electric 
field propagating along the z-axis) electromagnetic (EM) 
waves in a planar metamaterial or optical waveguide structure. 
In particular, we consider the case where the electric and mag- 
netic field components take the form E(x, z, t) = yE y (x, z, t) 
and H(x, z, t) = xH x (x, z, t) + zH z (x, z, t), where x, y, z are 
the unit vectors along the x, y, z directions, respectively, and 
we have assumed no variations (i.e., a homogeneous medium) 
with respect to the variable y. Under these assumptions, we 
may use Maxwell's equations - namely, Ampere's and Fara- 
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day's laws - take, respectively, the following form: 



dH z dH x _ dD y 
dx dz dt 



dE y 
~dz~ 



dz 
dt ' 



8Ey 

dx 



' dt 



(1) 

(2) 



Here D y is the v-component of the displacement vector D = 
yD y . Furthermore, we assume that the magnetic induction 
vector B is connected with the magnetic field intensity H by 
means of the constitutive relation B = fi(oj)H, where ft(aj) is 
the linear magnetic permeability (hereafter, we use / and / 
to denote any function / in the time- and frequency-domain, 
respectively). Additionally, we assume that the considered 
structure exhibits a weak cubic (Kerr-type) nonlinearity in 
its dielectric response. In other words, D y = e * E y + Pnl, 
where e is the permittivity, * denotes the convolution integral 
fit) * g(t) = £ + J f(r)g(t - r)dr of any functions f(t) and g(t), 
while the nonlinear polarization P NL is of the form, 

X+oo 
XNL(t-Ti,t-T 2 ,t-T 3 ) 
oo 

x Ey(Ti)Ey(T2)Ey(T 3 )dTidT 2 dT 3 . (3) 

Here, 6o is the dielectric constant of vacuum and xnl is the 
nonlinear electric susceptibility of the medium. In the case 
of small-amplitude, ultra-short pulse propagation, the non- 
linear response can safely be considered to be instantaneous, 
namely, 

XNL(t -T U t- T 2 , t ~ T 3 ) = K6(t ~ Ti)d(t - T 2 )S(t - T 3 ), (4) 

where k is the Kerr coefficient given by k = ±E~ 2 , with E c be- 
ing a characteristic electric field value; generally, both cases 
of focusing (k > 0) and defocusing (k < 0) dielectrics are pos- 
sible. Notice that Eqs. © and © imply that Pnl = £o*<E 3 and, 
thus, D y = e * E y + eoKEy. Substituting the considered form 
of the constitutive relations into Eqs. ©-©, we derive the 
following equation for the v-component of the electric field 
intensity (E y ) which, for convenience, will be denoted here- 
after by E\ 



V 2 E - d 2 (e *fi*E)- e Kd 2 (fi * E 3 ) = 0, 



(5) 



where V 2 = d 2 + d 2 is the Laplacian in the (x, z)-plane. 

Equation ([5]) is the (2 + l)-dimensional generalization of the 
ID Klein-Gordon t ype model derived in the context of nonlin- 
ear fiber optics 11121 12311 (in this case, /i = const.) and nonlinear 
metamaterials 113],|24]] (in this case, d^jl ± 0). Below, we will 
analyze the latter (more general) case, and assume that both 
permittivity and permeability are frequency dependent. In 
particular, considering general Kramers-Kronig or Sellmeier 
formulas (see, e.g., Ref. |22|]), we assume that the frequency 
dependence of e = e(co) and ft = fi(oj) can be approximated 
by the relations 



e(co) ~ 6 



(1(a)) ~ 



it-® 



(6) 



where a\, a 2 , J3\ and f5 2 are some constants. The above ap- 
proximations can be applied to the contexts of nonlinear left- 
handed metamaterials and nonlinear optical slab waveguides. 



Specifically, in the context of nonlinear left-handed metama- 
terials, a i = 1, a 2 = ci) 2 p ,Pi = l-F and ^2 = Fa) 2 es , where cj p , 
F and aj rQS denote, respectively, the plasma frequency, the fill- 
ing factor, and the magnetic permeability resonance frequency 
11311 . On the other hand, in the context of nonlinear optical 
slab waveguides, a\ = ef\ a 2 = ep\ J3\ = 1 and J3 2 = 0, 
where 6^ 0) and are relative dielectric constants (with 
being measured in units of squared angular frequency) ob- 
tained by matching the full form of the permittivity with the 
first of Eqs. © over a specific wavelength range in the infrared 
regime 111 211 . 

Next, we express Eq. © in the frequency domain and sub- 
stitute Eqs. ©, keeping terms up to order O(co~ 2 ) (i.e., assum- 
ing that ai/3 2 /co 4 « 1); then, expressing the resulting equation 
back in time domain, and measuring time, space, and field 
intensity E 2 in units of 1/ ^fa^, c/ ^a\a 2 f3\ and \k\~ 1 respec- 
tively, we reduce Eq. © in the following dimensionless form: 



V 2 E 



8 2 E 



aE ■ 



s K (J3E 3 +yd 2 E 3 ) 



0. 



(7) 



In the above equation, s K = sign(/<:), while the other constants 
are given by: 



1 ^ Pi 
a = — + . 

a\ a 2 pi 



a\a 2 P\ 



1 



(8) 



Note that in the context of nonlinear metamaterials a - 1 +fi = 
1 + Foj 2 es /[(l - F)qj 2 ] and y = 1, while in the context of 

nonlinear optical slab waveguides a = y= l/6^ 0) and/? = 0. 

We now consider propagation of ultrashort pulses, of width 
s, where < s «c 1 is a formal small parameter, which will 
also set the field amplitude (see below). Then, we employ the 
method of multiple scales to derive from Eq. © two different 
versions of short pulse equations in (2 + l)-dimensions. In that 
regard, we introduce the following asymptotic expansion for 
the unknown field E: 



E = sE\(T, X n , Z n ) + s 2 E 2 (T, X n , Z n ) + . . . , 



(9) 



where the functions E n depend on the spatial variables X n and 
Z n (n = 1, 2, . . .), as well as on the fast time variable T. Defin- 
ing Z n as: 



Z n = s n z, 



(10) 



we consider two different definitions for X n and T, namely: 

^ t-z 



(id 



£ 

where cr = ±1. Then, substituting Eqs. @, dTOj) and (1TTI) 
in Eq. ©, we derive at order 0(s) the following (2 + 1)- 
dimensional SPE for the unknown field E\ : 

d 2 E x d 2 Ei d 2 q 

2t^t t + <*E\ + s K y — AEl) = 0. (13) 

dZidT dX\ K ^0T 2 1 

Equation (TT3l) will be called hereafter SPE-I. Similarly, sub- 
stituting Eqs. ©, dT0l) and (fT2l) in Eq. ©, we derive [again at 
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order O(s)] another (2 + l)-dimensional version of the SPE, 
namely: 

^{^ + Sk) +aE ^^w^--^ (14) 

which will be called hereafter SPE-II. We note that variants of 
these models have been considered in the past in the context 
of ultrashort propagation in nonlinear dielectrics 021 |4|, and 
more recently in relevant studi es H , as well as in the context 
of collapse in two-level media lUlh . 

III. THE SPE-I: GENERAL PROPERTIES AND 
NUMERICAL STUDY 

In this section, we focus on the SPE-I which, for simplicity 
of notation, is expressed in the form: 

2E zt -E xx + aE + s(E 3 ) tt = 0, (15) 

where s = s K y and subscripts denote partial derivatives, with z 
being the evolution variable. Below, we will consider general 
properties of this equation and discuss its solutions. 

A. Properties and canonical structure 

First, we study the Hamiltonian structure of Eq. (TBI ). For 
this purpose, we integrate Eq. (1T3T) with respect to time t and, 
introducing the auxiliary field E = <p t , we express Eq. (TBI) as 
follows: 

20 tz -(p xx + a(f) + s(tf) t = 0. (16) 

Then, it can be verified that Eq. (TT61) can be obtained from the 
variational principle, with Lagrangian density: 

£ = -^(p 2 + S -(f>t + Mz-\(pl (17) 
From this Lagrangian density, we can derive the Hamiltonian: 

J-oc J-oo \0<Pz I 

as well as the momenta: 

X+oo p*+oo r\ p r*+oo p*+oo 

I -^(p t dtdx= I I tfdtdx, (19) 
oo U —co "tZ *J — OO sj —CO 

X+oo s>+oo o p p*+oo p»+oo 

—<f> x dtdx= MAtAx. (20) 

oo \J —co UCPz %J —OO %J —co 

Let us next consider the Fourier transform of Eq. (fT5l with 
respect to time t, which leads to the equation: 

(icj)E z = X -E xx -°^E- S -{icofE\ (21) 



The above equation implies that 

i IE — aE - \ 

E z (co u x,z) = — — + scdE 3 , for oj ± 0, (22) 

2\ co ) 

E xx -aE = 0, for co = 0. (23) 

In our numerical simulations below, we will set E(co = 
0, x, z) - 0, so that Eq. (l23t is satisfied. Thus, in this case, the 
Fourier transform of the field E leads to the following "zero- 
mass constraint": 




E(t, x, z)dt = for any x, z, (24) 



which also holds for the traditional SPE in (1 + l)-dimensions 
(see, e.g., the relevant analysis of Ref. i27ll ). 

B. ID breather-like structures 

Let us now seek one-dimensional (ID) solutions of 
Eq. (TBI) , by assuming that the unknown field E depends on 
the traveling- wave coordinates £ and 77, defined as: 

c 2 

f = z, T] = t + cx+ —z, (25) 

where c is an arbitrary real constant setting the velocities of 
the traveling wave in the (x, t) and (x, z) planes. Using the 
above variables, Eq. (TT3T) is reduced to the form: 

IE & + aE + s(E 3 ) VT1 = 0, (26) 

which is actually the (1 + l)-dimensional SPE fl3l . As shown 
in the simulations of Ref. 111811 . the most robust among the 
various solutions of the ID SPE is the breather-like struc- 
ture (this solution satisfies the zero-mass constraint). Nat- 
urally, this purely ID structure satisfies the full 2D SPE-I, 
Eq. (IT5l> and, thus, an interesting question concerns the sta- 
bility of this solution in the 2D space. A similar question ap- 
pears in many physically relevant models and, in many cases, 
the answer is that such "planar" solutions are prone to trans- 
verse instabilities in higher-dimensional settings; as charac- 
teristic examples, we mention the line soliton solutions of the 
Kadomtsev-Petviashvilli-I (KP-I) equation which decay into 
lumps |28l|29|], or the dark soliton stripes of the defocusing 
NLS equation which decay into vortices (3(1 HH in (2 + 1)- 
dimensions. 

To study the stability of the ID breather-like solution of 
Eq. (l26l) in the framework of the the full 2D SPE-I, we have 
used the following procedure. We employed the breather so- 
lution of the ID SPE Eq. (f26l) and also added, as a perturba- 
tion, a small noisy signal, of amplitude of 1% of the breather's 
amplitude. Then, using the resulting structure (cf. top panel 
of Fig. [T]) as an initial condition, we numerically integrated 
Eq. (1131) by means of a Galerkin method (and assuming pe- 
riodic boundary conditions in our numerical scheme). The 
results (corresponding to parameter values c = 0, s = -1/3, 
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FIG. 1: (Color online) Top (four) panels: Contour plots showing 
profiles of the evolution of a perturbed ID breather in the (x, t) plane, 
when evolved according to Eq. (TBI ). Snapshots correspond to z = 
(top left), z = 50 (top right) z = 100 (bottom left), and z = 150 
(bottom right). Bottom panel: The evolution of the breather for x = 
0. Parameter values are c = 0, s = -1/3, and a = -2. 
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FIG. 2: (Color online) The evolution of the conserved quantities of 
the SPE-I for the simulation shown in Fig. [T] Top panel shows the 
Hamiltonian H and bottom panels show the momenta M t (left) and 
M t (right). The relative error for H and M t is of order 10~ 4 and 10~ 3 
respectively, and M x is zero. 



and a = -2) are shown in the panels of Fig. [T] in terms of 
different contour plots depicting the profile of the ID breather 
in the (x, t) plane for various values of the propagation dis- 
tance z (and also the evolution at x = as a function of (z, t)). 
It is clear that the breather is robust, at least up to z = 200 
(where the simulation ended). We should also mention that for 
these simulations, we have also calculated the evolution of the 
Hamiltonian and momenta [cf. Eqs. ([T8l)-([20bl. The results, 
depicted in Fig. [2 justify the conservation of these quantities 
with a satisfactory (relative) accuracy, of order 10" 3 or less. 

We finally note that similar results (not shown here) were 
also obtained for breathers with i.e., for c ± in Eq. (|25l) ; in 
such a case, the only difference is that the breather is "tilted", 
i.e., oblique in the (x, t) plane with respect to its direction in 
the case c = and follows a similar evolution (i.e., it is stable 
up to end of the simulation time). 

C. Localized initial data 

Having discussed the properties of the ID breather in the 
2D setting, we now turn our attention to initial data asso- 
ciated with Eq. ([26]), which are localized in both transverse 
directions, x and t. In that regard, it is convenient to con- 
sider at first the decomposition E(x, t) = f(t)g(x), and sub- 
stitute this ansatz in the zero-mass constraint, Eq. (l24t . This 
way, for nontrivial solutions, we derive the necessary condi- 
tion JT/(*)A = 0. 

Taking into account the above constraint, we now may use 
fit) = (1 - t 2 )exp(-t 2 /2) which has the above property, and 
also choose g(x) to be of the same functional form, namely: 
g(x) = (1 - x 2 )exp(-x 2 /2). The aim of the latter choice is 
to produce a two-dimensional localized waveform. Employ- 
ing these choices, we can now numerically integrate Eq. ([T5l . 
using the initial condition: 

E(z = 0, x, t) = (1 - x 2 )(l - t 2 ) exp[-(x 2 + t 2 )/2]. (27) 

The results of our simulations are presented in Fig. [3] where 
we show the evolution of this initial data. It is clearly ob- 
served that, already at small values of the propagation dis- 
tance (z ~ 2), the initially localized structure bends and splits 
at (x, t) = (0, 0), thus forming two "wing-like" structures. The 
size (length) of these structures is small at the early stages of 
the evolution but, afterwards, their spatial extent is increased, 
as the initial data progressively disperses. This way, the result- 
ing structures yield, at longer propagation distances (see bot- 
tom panels of Fig. [3}, a quasi-one-dimensional pattern, some- 
what reminiscent of the breather states examined previously. 
Here we should mention that our simulations end up at rela- 
tively small distance (z = 38) in order to avoid interference of 
these expanding quasi- ID structures with the boundaries (re- 
call that we use periodic boundary conditions in our numerical 
scheme). 

We also note in passing that we have tried other localized 
initial conditions, which led to qualitatively similar results: in 
all cases, the respective evolutions of the initial localized data 
gradually transformed into quasi-one-dimensional dispersing 
structures of the above type. 
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FIG. 3: (Color online) Contour plots showing profiles of the field 
E, in the {x, t) plane, evolved as per SPE-I, Eq. (TT5K with localized 
initial data [cf. Eq. (ETVl. The snapshots, from left to right, and top 
to bottom correspond to z = 0, 3, 10, 20, 30, 38. Parameter values are 
.v = -1/3 and a = -2. 



IV. THE SPE-II: GENERAL PROPERTIES AND 
NUMERICAL STUDY 

Let us now consider the SPE-II which we express, for sim- 
plicity of notation, in the following form: 



2a(E zt + E xt ) + aE + s(E 3 ) tt = 0, 



(28) 



where s = s K y, as in the SPE-I. Below we will follow 
the presentation of the previous section and discuss gen- 
eral properties of this model, such as the corresponding 
Lagrangian/Hamiltonian structure and relevant conservation 
laws as well as some of its prototypical solutions. 



A. Properties and canonical structure 

As in the case of the SPE-I, we integrate Eq. (l28t with re- 
spect to t and, introducing the field E = cp t , we express SPE-II 
in the following form: 



2cr((p xt + cp zt ) + acp + s((f> 3 t ) t = 0. 



(29) 



The above equation can be obtained from the variational prin- 
ciple, with Lagrangian density: 

£ = ~\ 02 + + + Mx) ' (30) 



The corresponding Hamiltonian can then be found as: 

H - £X(IM d,d * 

X+oo ~+oo 
J {2<l> 2 - 4$ - o-Mxjdtdx, (31) 



while the momenta read: 



I - — (p t dtdx = cr I I tfdtdx, 

oo sJ — oo ^YZ %J—oo U -co 

(32) 

X+oo a»+oo Q£ r*+oo r*+co 

\ -WT-<f>xdtdx = cr I I (p t (p x dtdx. 
oo *J — oo ^tZ, %J —OO xJ —oo 

(33) 

Next, we consider the Fourier transform of Eq. (|28l with 
respect to t, which leads to the equation: 

2o-(ioj)(E z + E x ) + aE + s(E 3 )(ico) 2 = 0. (34) 
Solving the above equation with respect to E we find: 



E z (a), x, z) = -E x - 



s(ico) 

2cr(ioj) 2cr 



-E - 



(£ 3 ), for co ^ 0, (35) 



E(oj = 0,x,z) = 0, fora) ± 0. 



(36) 



The latter equation leads again to the zero-mass constraint [cf . 
Eq. (l24l) l that we found in the case of the SPE-I as well. This 
condition will also be satisfied in our simulations below. 



B. ID breathers and initial data localized in 2D 

We consider traveling wave solutions of Eq. (I28K in the 
form £(£, rj), where the coordinates £ and r] are defined as: 



t + cx - cz, 



(37) 



where c is an arbitrary real constant. This way, Eq. 
transformed to the equation: 



is 



IcrEfr + aE + s(E \ 



0, 



(38) 



which is actually the ID SPE model Ill2h . Since the latter 
admits breather solutions, we may follow the procedure de- 
scribed in the previous section and study numerically the evo- 
lution of such a ID solution in the 2D setting of Eq. d28l) . 

In this case also, the numerical integration of Eq. (l28t has 
shown that this ID solution is stable in the 2D setting (as was 
also in the framework of the SPE-I). An example (pertaining 
to parameter values c = 0, s = -1/3, a = -2, and cr = 1) 
is shown in Fig.|U Additionally, the numerical calculation of 
the evolution of the Hamiltonian and momenta [cf. Eqs. (f3TT> - 
(l33t l depicted in Fig. [5] illustrate the conservation of these 
quantities with a relative error of order 10" 4 or less. We also 
note that similar results (not shown here) were also obtained 
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FIG. 5: (Color online) The evolution of the conserved quantities of 
SPE-II for the simulation shown in Fig. U) the top panel shows the 
Hamiltonian H and the bottom panels show the momenta M t (left) 
and M x (right). The relative error for H and M t is of order 10~ 4 
or less, and the value of M x is zero (and remains so throughout the 
simulation). 



FIG. 4: (Color online) Top (four) panels: Contour plots showing pro- 
files of the ID breather solution of Eq. (l38l in the (x, t) plane, when 
evolved as per Eq. (|28] >. Snapshots correspond to z = (top left), 
z = 10 (top right) z = 30 (bottom left), and z = 50 (bottom right). 
Bottom panel: The evolution of the breather for x = 0. Parameter 
values are: c = 0, s = -1/3, a = -2, and cr - 1 . 



for oblique moving breathers, i.e., for c ^ in Eq. dTTb , as in 
the case of SPE-I. 

Finally, as in the case of SPE-I, we study the evolution of 
localized data in 2D (i.e., in both x and t) in the framework 
of the SPE-II. A typical example of the result obtained by 
the numerical integration of Eq. (TBI ) with such localized ini- 
tial data is shown in Fig. [6| (parameter values are s = -1/3, 
a = -2, and cr = 1). It is observed that after a small prop- 
agation distance (z « 5) the initially localized waveform be- 
gins to broaden along the r-axis, but still remaining localized 
along the x-axis. As a result, a quasi- ID structure is gradually 
formed, which travels faster along the x-direction (where it is 
localized) than in the t direction (where it is elongated). In 
the latter direction, the structure also possesses an alternating 
spatial structure which merits further investigation. 



V. DISCUSSION AND CONCLUSIONS 

In conclusion, we have derived from Maxwell's equations 
two (2 + l)-dimensional short pulse equations, referred to as 
SPE-I and SPE-II. These equations may find applications in 
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FIG. 6: (Color online) Contour plots showing profiles of the field 
E, in the (x, t) plane, evolved as per SPE-II, Eq. dTSb , with localized 
initial data [cf. Eq. (12711. The snapshots, from top to bottom and 
left to right, correspond to z = 0,3, 10,20,30,40. The domain for 
numerical computation in xx t plane is [-20tt, 207r] X [-20tt, 20tt], 
we are zooming in here in the snapshots to show more detail. 
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various physical contexts where the study of ultrashort elec- 
tromagnetic pulses is important; such contexts include non- 
linear metamaterials, nonlinear optical waveguide structures, 
nonlinear dielectric media, and others. Since both SPE-I 
and SPE-II actually generalize the (1 + l)-dimensional SPE 
fl2h . they can be used to the study of transverse (diffraction- 
induced) dynamics of ultra-short pulses in such settings. Suit- 
able assumptions on the nature of the electric and magnetic 
field and the form of the permittivity and permeability under 
which the equations can be derived were provided. 

We have found and presented various general properties of 
SPE-I and SPE-II. Particularly, we have identified the La- 
grangian and Hamiltonian structure, and have used invari- 
ances to infer (from Noether's theory) the corresponding mo- 
menta, as well as the associated zero-mass constraints; the 
latter, have to be satisfied for the solutions of these equations 
and, thus, are also associated with the choice of the initial 
data used for the numerical integration of SPE-I and SPE-II. 
We have conducted a series of numerical experiments for the 
2D SPEs using, as initial conditions, either the ID breather 
solution of the underlying ID SPEs or a localized (in 2D) 
waveform - both satisfying the zero-mass constraint. Our mo- 
tivation was to study the stability and transverse dynamics of 
the most robust solutions of the ID analogue of the system, 
and also examine the fate of purely 2D initial data and poten- 
tially identify structures that can be supported by the SPE-I 
and SPE-II models. 

Our numerical simulations have shown that the ID 
breathers propagate (even when they are initially perturbed 
by a small noise) practically undistorted. An important con- 
clusion is that these ultrashort localized structures are actually 
insensitive in the presence of diffraction or, in other words, 
they appear to be robust in the presence of (small) transverse 



perturbations for propagation distances of the order of a few 
hundred dimensionless units. On the other hand, simulations 
employing initial data localized in 2D have shown that, during 
evolution, the initial data gradually transforms into quasi- ID 
structures (which differ between SPE-I and SPE-II). In fact, 
we were not able to find any, purely 2D, nonlinear waveform 
that can be supported by either the SPE-I or the SPE-II. 

The above results were obtained in the framework of the 
particular models, i.e., SPE-I and SPE-II, that we derived 
and considered in this work. It would be interesting to per- 
form similar studies (i.e., transverse dynamics of ID ultrashort 
pulses and localized 2D structures) in the framework of other 
versions of the SPE-I and SPE-II, stemming from the incorpo- 
ration of higher-order effects (as in the ID case, in the context 
of the so-called regularized SPE J23l|32|]). On the other hand, 
still in the context of SPE-I and SPE-II, it would be relevant 
to consider other types of solutions, e.g., loop-type solutions 
or periodic waveforms composed by breathers or loom (as in 
the spirit of the analysis in the ID case - see Ref. [18]), and 
others. It would also be relevant to compare the properties of 
the models derived herein with those of other models for ul- 
trashort pulses including, e.g., iBl- lUlll . Finally, generalizing 
the present models to (3 + l)-dimensions, removing the as- 
sumption of spatial homogeneity along the y-direction would 
also constitute an interesting theme for future studies. 
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